Ejercicio de Simulación
Vamos a hacer un ejercicio de simulación para ver como se identifica
la componente estacional.
library(urca)
library(forecast)
library(tseries)
library(lmtest)
library(uroot)
library(fUnitRoots)
Loading required package: timeDate
Loading required package: timeSeries
Attaching package: ‘timeSeries’
The following object is masked from ‘package:zoo’:
time<-
Loading required package: fBasics
Attaching package: ‘fUnitRoots’
The following objects are masked from ‘package:urca’:
punitroot, qunitroot, unitrootTable
library(sarima)
Loading required package: stats4
Attaching package: ‘sarima’
The following object is masked from ‘package:stats’:
spectrum
require("PolynomF")
Loading required package: PolynomF
###Simulación de un proceso con raíz unitaria estacional
#x11()
x <- ts(sarima::sim_sarima(n=144, model = list(iorder=0, siorder=1, nseasons=12, sigma2 = 1),n.start=24),frequency = 12)
plot(x)

acf(x,lag.max = 36)

monthplot(x)

nsdiffs(x)####Decreta cuantas diferencias estacional a través de la aplicación de
[1] 1
###Algunas pruebas de raíces unitarias estacionales.
###diferencia estacional
Dx=diff(x,lag=12,differences = 1)###lag:periodo s.
plot(Dx)

acf(Dx,lag.max = 36)

monthplot(Dx)

nsdiffs(Dx)
[1] 0
####Simulación de un SAR
#x11()
x1 <- ts(sim_sarima(n=144, model = list(ar=c(rep(0,11),0.8)),n.start=24),frequency=12)
plot(x1)

acf(x1,lag.max = 36)

monthplot(x1)

nsdiffs(x1)
[1] 1
ndiffs(x1)
[1] 1
p <- polynom(c(1,c(rep(0,11),-0.8)))
solve(p)
[1] -1.0187693+0.0000000i -0.8822801-0.5093846i
[3] -0.8822801+0.5093846i -0.5093846-0.8822801i
[5] -0.5093846+0.8822801i 0.0000000-1.0187693i
[7] 0.0000000+1.0187693i 0.5093846-0.8822801i
[9] 0.5093846+0.8822801i 0.8822801-0.5093846i
[11] 0.8822801+0.5093846i 1.0187693+0.0000000i
abs(solve(p))
[1] 1.018769 1.018769 1.018769 1.018769 1.018769 1.018769
[7] 1.018769 1.018769 1.018769 1.018769 1.018769 1.018769
###Note lo cerca que están la raíces de la no estacionariedad del proceso, por eso
####aunque si bien el proceso es estacionario, notamos hay una cercanía a
####e tener una compoenete estacional.
####El anterior modelo puede escribirse como:
x2 <- ts(sim_sarima(n=144, model=list(sar=0.8, iorder=0, siorder=0, nseasons=12),n.start=24),frequency = 12)
plot(x2)

acf(x2, lag.max=48)

monthplot(x2)

nsdiffs(x2)
[1] 1
Ejemplo Pasajeros
Vamos a ver como se hace el modelamiento completo de la serie de
pasajeros.
Iniciaremos con la transformación Box-Cox y las pruebas de raíces
Unitarias


[1] 1
Call:
ar(x = lAirpassengers)
Coefficients:
1 2 3 4 5 6
1.0013 -0.0959 0.0329 -0.0252 0.0231 0.0005
7 8 9 10 11 12
-0.0157 -0.0762 0.1071 -0.0236 0.0675 0.4536
13
-0.4854
Order selected 13 sigma^2 estimated as 0.01336
Augmented Dickey-Fuller Test
data: lAirpassengers
Dickey-Fuller = -2.147, Lag order = 13, p-value =
0.5152
alternative hypothesis: stationary
Warning in fUnitRoots::adfTest(lAirpassengers, lags = 12, type = "nc") :
p-value greater than printed p-value
Title:
Augmented Dickey-Fuller Test
Test Results:
PARAMETER:
Lag Order: 12
STATISTIC:
Dickey-Fuller: 3.7872
P VALUE:
0.99
Description:
Thu Jun 9 11:01:18 2022 by user:
###############################################
# Augmented Dickey-Fuller Test Unit Root Test #
###############################################
Test regression none
Call:
lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
Residuals:
Min 1Q Median 3Q Max
-0.103864 -0.023652 -0.001455 0.022160 0.126649
Coefficients:
Estimate Std. Error t value Pr(>|t|)
z.lag.1 0.005439 0.001436 3.787 0.000241 ***
z.diff.lag1 -0.198813 0.072080 -2.758 0.006737 **
z.diff.lag2 -0.273730 0.073309 -3.734 0.000292 ***
z.diff.lag3 -0.233607 0.072258 -3.233 0.001589 **
z.diff.lag4 -0.293133 0.073926 -3.965 0.000126 ***
z.diff.lag5 -0.206562 0.072058 -2.867 0.004915 **
z.diff.lag6 -0.266919 0.071493 -3.734 0.000292 ***
z.diff.lag7 -0.234526 0.071847 -3.264 0.001436 **
z.diff.lag8 -0.327393 0.073197 -4.473 1.79e-05 ***
z.diff.lag9 -0.198455 0.073623 -2.696 0.008054 **
z.diff.lag10 -0.279931 0.072710 -3.850 0.000192 ***
z.diff.lag11 -0.176122 0.073011 -2.412 0.017394 *
z.diff.lag12 0.627402 0.072683 8.632 3.34e-14 ***
---
Signif. codes:
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.04308 on 118 degrees of freedom
Multiple R-squared: 0.8564, Adjusted R-squared: 0.8406
F-statistic: 54.14 on 13 and 118 DF, p-value: < 2.2e-16
Value of test-statistic is: 3.7872
Critical values for test statistics:
1pct 5pct 10pct
tau1 -2.58 -1.95 -1.62
Title:
Augmented Dickey-Fuller Test
Test Results:
PARAMETER:
Lag Order: 12
STATISTIC:
Dickey-Fuller: -1.5325
P VALUE:
0.7711
Description:
Thu Jun 9 11:01:18 2022 by user:


Call:
ar(x = dlAirpassengers)
Coefficients:
1 2 3 4 5 6
-0.0993 -0.1226 -0.2824 -0.2875 -0.1184 -0.2406
7 8 9 10 11 12
-0.1829 -0.2971 -0.1040 -0.2627 -0.1259 0.5734
13 14 15
0.0179 -0.1667 0.1200
Order selected 15 sigma^2 estimated as 0.002988
Warning in tseries::adf.test(dlAirpassengers, k = 2) :
p-value smaller than printed p-value
Augmented Dickey-Fuller Test
data: dlAirpassengers
Dickey-Fuller = -7.6877, Lag order = 2, p-value =
0.01
alternative hypothesis: stationary
Title:
Augmented Dickey-Fuller Test
Test Results:
PARAMETER:
Lag Order: 12
STATISTIC:
Dickey-Fuller: -1.2228
P VALUE:
0.2243
Description:
Thu Jun 9 11:01:18 2022 by user:
Identificación de la componente ARMA estacional y la componente ARMA
ordinaria
####################################
library(uroot)
require(forecast)
######Diferencia Estacional(continuación AirPassengers)#######
monthplot(dlAirpassengers)
nsdiffs(dlAirpassengers)
[1] 1
nsdiffs(AirPassengers)
[1] 1
DdlAirpassengers=diff(dlAirpassengers,lag=12)###lag=s
#x11()
par(mfrow=c(2,1))

plot(dlAirpassengers)
plot(DdlAirpassengers)

monthplot(DdlAirpassengers)
nsdiffs(DdlAirpassengers)
[1] 0
##Autocorrelogramas
#x11()
acf(DdlAirpassengers)

acf(DdlAirpassengers,lag.max = 48, ci.type='ma')# q=0,1, Q=0,1
pacf(DdlAirpassengers,lag.max = 48) # p=0,1,2,...,9, P=0,1

#SARIMA(p=0,d=1,q=1)x(P=0,D=1,Q=1)s=12
Ajuste del Modelo y Análisis de Residuales
##Ajuste del modelo
###Arima Estacional o SARIMA(p=0,d=1,q=3)x(P=0,D=1,Q=1)s=12 con transformación logaritmica
#Modelo MA(1) estacional
modelo = Arima(AirPassengers, c(0, 1, 1),seasonal = list(order = c(0, 1, 1), period = 12),lambda = 0)
coeftest(modelo)
z test of coefficients:
Estimate Std. Error z value Pr(>|z|)
ma1 -0.401828 0.089644 -4.4825 7.378e-06 ***
sma1 -0.556945 0.073100 -7.6190 2.557e-14 ***
---
Signif. codes:
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
modeloalter= Arima(AirPassengers, c(1, 1, 0),seasonal = list(order = c(1, 1, 0), period = 12),lambda = 0)
## Análisis de residuales
#x11()
residuales <- modelo$residuals
plot(residuales)

acf(residuales,lag.max = 24)

pacf(residuales,lag.max = 24)
#Test de autocorrelaci?n
Box.test(residuales, lag = (length(residuales)/4), type = "Ljung-Box", fitdf = 2)
Box-Ljung test
data: residuales
X-squared = 37.874, df = 34, p-value = 0.2969
######Análisis de Outliers
#Test de normalidad
jarque.bera.test(residuales)
Jarque Bera Test
data: residuales
X-squared = 5.2265, df = 2, p-value = 0.0733
###Estad?ticas CUSUM
res=residuales
cum=cumsum(res)/sd(res)
N=length(res)
cumq=cumsum(res^2)/sum(res^2)
Af=0.948 ###Cuantil del 95% para la estad?stica cusum
co=0.14422####Valor del cuantil aproximado para cusumsq para n/2
LS=Af*sqrt(N)+2*Af*c(1:length(res))/sqrt(N)
LI=-LS
LQS=co+(1:length(res))/N
LQI=-co+(1:length(res))/N
par(mfrow=c(2,1))

plot(cum,type="l",ylim=c(min(LI),max(LS)),xlab="t",ylab="",main="CUSUM")
lines(LS,type="S",col="red")
lines(LI,type="S",col="red")
#CUSUM Square
plot(cumq,type="l",xlab="t",ylab="",main="CUSUMSQ")
lines(LQS,type="S",col="red")
lines(LQI,type="S",col="red")

Pronóstico
ecm
[1] 342.0955

Vale la pena decir que si los residuales son no correlacionados pero
ellos no tiene distribución normal, se puede usar los intervalos
bootatrap ver sección 5.5 libro fpp3. Esto se hace añadiendo
bootstrap=TRUE en la función forecast().
LS0tCnRpdGxlOiAiU0FSSU1BIgpvdXRwdXQ6IAogIGdpdGh1Yl9kb2N1bWVudDogZGVmYXVsdAogIGh0bWxfbm90ZWJvb2s6IGRlZmF1bHQKLS0tCgpgYGB7ciBzZXR1cCwgaW5jbHVkZT1GQUxTRX0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KGVjaG8gPSBUUlVFKQpgYGAKCiMjIEVqZXJjaWNpbyBkZSBTaW11bGFjacOzbgoKVmFtb3MgYSBoYWNlciB1biBlamVyY2ljaW8gZGUgc2ltdWxhY2nDs24gcGFyYSB2ZXIgY29tbyBzZSBpZGVudGlmaWNhIGxhIGNvbXBvbmVudGUgZXN0YWNpb25hbC4KCmBgYHtyIFNpbXVsYWNpb259CmxpYnJhcnkodXJjYSkKbGlicmFyeShmb3JlY2FzdCkKbGlicmFyeSh0c2VyaWVzKQpsaWJyYXJ5KGxtdGVzdCkKbGlicmFyeSh1cm9vdCkKbGlicmFyeShmVW5pdFJvb3RzKQpsaWJyYXJ5KHNhcmltYSkKcmVxdWlyZSgiUG9seW5vbUYiKQoKIyMjU2ltdWxhY2nDs24gZGUgdW4gcHJvY2VzbyBjb24gcmHDrXogdW5pdGFyaWEgZXN0YWNpb25hbAojeDExKCkKeCA8LSB0cyhzYXJpbWE6OnNpbV9zYXJpbWEobj0xNDQsIG1vZGVsID0gbGlzdChpb3JkZXI9MCwgc2lvcmRlcj0xLCBuc2Vhc29ucz0xMiwgc2lnbWEyID0gMSksbi5zdGFydD0yNCksZnJlcXVlbmN5ID0gMTIpCnBsb3QoeCkKYWNmKHgsbGFnLm1heCA9IDM2KQptb250aHBsb3QoeCkKbnNkaWZmcyh4KSMjIyNEZWNyZXRhIGN1YW50YXMgZGlmZXJlbmNpYXMgZXN0YWNpb25hbCBhIHRyYXbDqXMgZGUgbGEgYXBsaWNhY2nDs24gZGUgCiMjI0FsZ3VuYXMgcHJ1ZWJhcyBkZSByYcOtY2VzIHVuaXRhcmlhcyBlc3RhY2lvbmFsZXMuCgoKIyMjZGlmZXJlbmNpYSBlc3RhY2lvbmFsCkR4PWRpZmYoeCxsYWc9MTIsZGlmZmVyZW5jZXMgPSAxKSMjI2xhZzpwZXJpb2RvIHMuCnBsb3QoRHgpCmFjZihEeCxsYWcubWF4ID0gMzYpCm1vbnRocGxvdChEeCkKbnNkaWZmcyhEeCkKIyMjI1NpbXVsYWNpw7NuIGRlIHVuIFNBUgojeDExKCkKeDEgPC0gdHMoc2ltX3NhcmltYShuPTE0NCwgbW9kZWwgPSBsaXN0KGFyPWMocmVwKDAsMTEpLDAuOCkpLG4uc3RhcnQ9MjQpLGZyZXF1ZW5jeT0xMikKcGxvdCh4MSkKYWNmKHgxLGxhZy5tYXggPSAzNikKbW9udGhwbG90KHgxKQpuc2RpZmZzKHgxKQpuZGlmZnMoeDEpCnAgPC0gcG9seW5vbShjKDEsYyhyZXAoMCwxMSksLTAuOCkpKQpzb2x2ZShwKQphYnMoc29sdmUocCkpCiMjI05vdGUgbG8gY2VyY2EgcXVlIGVzdMOhbiBsYSByYcOtY2VzIGRlIGxhIG5vIGVzdGFjaW9uYXJpZWRhZCBkZWwgcHJvY2VzbywgcG9yIGVzbwojIyMjYXVucXVlIHNpIGJpZW4gZWwgcHJvY2VzbyBlcyBlc3RhY2lvbmFyaW8sIG5vdGFtb3MgaGF5IHVuYSBjZXJjYW7DrWEgYSAKIyMjI2UgdGVuZXIgdW5hIGNvbXBvZW5ldGUgZXN0YWNpb25hbC4KIyMjI0VsIGFudGVyaW9yIG1vZGVsbyBwdWVkZSBlc2NyaWJpcnNlIGNvbW86CngyIDwtIHRzKHNpbV9zYXJpbWEobj0xNDQsIG1vZGVsPWxpc3Qoc2FyPTAuOCwgaW9yZGVyPTAsIHNpb3JkZXI9MCwgbnNlYXNvbnM9MTIpLG4uc3RhcnQ9MjQpLGZyZXF1ZW5jeSA9IDEyKQpwbG90KHgyKQphY2YoeDIsIGxhZy5tYXg9NDgpCm1vbnRocGxvdCh4MikKbnNkaWZmcyh4MikKYGBgCgojIyBFamVtcGxvIFBhc2FqZXJvcwoKVmFtb3MgYSB2ZXIgY29tbyBzZSBoYWNlIGVsIG1vZGVsYW1pZW50byBjb21wbGV0byBkZSBsYSBzZXJpZSBkZSBwYXNhamVyb3MuCgpJbmljaWFyZW1vcyBjb24gbGEgdHJhbnNmb3JtYWNpw7NuIEJveC1Db3ggeSBsYXMgcHJ1ZWJhcyBkZSByYcOtY2VzIFVuaXRhcmlhcwoKYGBge3IgcGFzc2VuZ2VycywgZWNobz1GQUxTRX0KIyMjIyMjQWp1c3RlIFNlcmllIGRhdG9zIEFpclBhc3NlbmdlcnMgcG9yIG1lZGlvIGRlbAojbW9kZWxvIEFSSU1BIGVzdGFjaW9uYWwsIGFzw60gY29tbyBzdSBjb3JyZXNwb25kaWVudGUKI2Fuw6FsaXNpcyBkZSByZXNpZHVhbGVzIHkgcHJvbsOzc3RpY29zCgpkYXRhKEFpclBhc3NlbmdlcnMpCnBsb3QoQWlyUGFzc2VuZ2VycykKIyMjIyMjIyMjI0Rlc3B1w6lzIGRlIGhhYmVyIGFwbGljYWRvIGxhIGRpZmVyZW5jaWEgb3JkaW5hcmlhIHkgZXN0YWNpb25hbAojIyMjIyNQcm9jZWRlbW9zIGEgdHJhdGFyIGRlIGlkZW50aWZpY2FyIGxhIGVzdHJ1Y3R1cmEgZGUgYXV0b2NvcnJlbGFjacOzbgojIyNhIGNvcnRvIHBsYXpvKEFSTUEpIHkgZXN0YWNpb25hbCBTQVJNQS4KIyMjUGFyYSBlc28sIGVzIG5lY2VzYXJpbyBoYWJlciBjb252ZXJ0aWRvIGxhIHNlcmllIGEgZXN0YWNpb25hcmlhIyMjCmxBaXJwYXNzZW5nZXJzPWxvZyhBaXJQYXNzZW5nZXJzKQpwbG90KGxBaXJwYXNzZW5nZXJzKQptb250aHBsb3QobEFpcnBhc3NlbmdlcnMpCgpmb3JlY2FzdDo6bmRpZmZzKGxBaXJwYXNzZW5nZXJzKQojIyMjI1BydWViYSBkZSBEaWNrZXkgRnVsbGVyIyMjIyMjCmFyKGxBaXJwYXNzZW5nZXJzKQp0c2VyaWVzOjphZGYudGVzdChsQWlycGFzc2VuZ2VycyxrPTEzKQpmVW5pdFJvb3RzOjphZGZUZXN0KGxBaXJwYXNzZW5nZXJzLGxhZ3MgPSAxMix0eXBlPSduYycpICAgIyMjSGF5IGxhIHByZXNlbmNpYSBkZSBSYcOteiBVbml0YXJpYQpzdW1tYXJ5KHVyY2E6OnVyLmRmKGxBaXJwYXNzZW5nZXJzLCBsYWdzID0gMTIpKQphZGZUZXN0KGxBaXJwYXNzZW5nZXJzLGxhZ3M9MTIsdHlwZT0nY3QnKSAgIyMjI1B1ZWRlIHRhbWJpw6luIGluZGljYXIKIyMjI0xhIHByZXNlbmNpYSBkZSB1bmEgdGVuZGVuY2lhIGRldGVybWluw61zdGljYQoKIyMjI0RpZmVyZW5jaWEgT3JkaW5hcmlhIyMjIyMjIyMjIyMjCgpkbEFpcnBhc3NlbmdlcnM9ZGlmZihsQWlycGFzc2VuZ2VycyxsYWc9MSkKI3gxMSgpCnBhcihtZnJvdz1jKDIsMSkpCnBsb3QobEFpcnBhc3NlbmdlcnMpCnBsb3QoZGxBaXJwYXNzZW5nZXJzKQphcihkbEFpcnBhc3NlbmdlcnMpCgp0c2VyaWVzOjphZGYudGVzdChkbEFpcnBhc3NlbmdlcnMsayA9IDIpICAjIyNObyBzZSBkZWJlIGRpZmVyZW5jaWFyIG3DoXMjIyMKCmZVbml0Um9vdHM6OmFkZlRlc3QoZGxBaXJwYXNzZW5nZXJzLGxhZ3MgPSAxMix0eXBlPSduYycpIAoKYGBgCgojIyBJZGVudGlmaWNhY2nDs24gZGUgbGEgY29tcG9uZW50ZSBBUk1BIGVzdGFjaW9uYWwgeSBsYSBjb21wb25lbnRlIEFSTUEgb3JkaW5hcmlhCgpgYGB7ciBDb21wb25lbnRlIEVzdGFjaW9uYWwgeSBPcmRpbmFyaWF9CiMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIyMjIwpsaWJyYXJ5KHVyb290KQpyZXF1aXJlKGZvcmVjYXN0KQojIyMjIyNEaWZlcmVuY2lhIEVzdGFjaW9uYWwoY29udGludWFjacOzbiBBaXJQYXNzZW5nZXJzKSMjIyMjIyMKbW9udGhwbG90KGRsQWlycGFzc2VuZ2VycykKbnNkaWZmcyhkbEFpcnBhc3NlbmdlcnMpCm5zZGlmZnMoQWlyUGFzc2VuZ2VycykKCkRkbEFpcnBhc3NlbmdlcnM9ZGlmZihkbEFpcnBhc3NlbmdlcnMsbGFnPTEyKSMjI2xhZz1zCiN4MTEoKQpwYXIobWZyb3c9YygyLDEpKQpwbG90KGRsQWlycGFzc2VuZ2VycykKcGxvdChEZGxBaXJwYXNzZW5nZXJzKQptb250aHBsb3QoRGRsQWlycGFzc2VuZ2VycykKbnNkaWZmcyhEZGxBaXJwYXNzZW5nZXJzKQoKCiMjQXV0b2NvcnJlbG9ncmFtYXMKI3gxMSgpCmFjZihEZGxBaXJwYXNzZW5nZXJzKQphY2YoRGRsQWlycGFzc2VuZ2VycyxsYWcubWF4ID0gNDgsIGNpLnR5cGU9J21hJykjIHE9MCwxLCBRPTAsMQpwYWNmKERkbEFpcnBhc3NlbmdlcnMsbGFnLm1heCA9IDQ4KSAjIHA9MCwxLDIsLi4uLDksIFA9MCwxCiNTQVJJTUEocD0wLGQ9MSxxPTEpeChQPTAsRD0xLFE9MSlzPTEyCmBgYAoKCiMjIEFqdXN0ZSBkZWwgTW9kZWxvIHkgQW7DoWxpc2lzIGRlIFJlc2lkdWFsZXMKYGBge3IgQWp1c3RlcyB5IFJlc2lkdWFsZXN9CiMjQWp1c3RlIGRlbCBtb2RlbG8KIyMjQXJpbWEgRXN0YWNpb25hbCBvIFNBUklNQShwPTAsZD0xLHE9Myl4KFA9MCxEPTEsUT0xKXM9MTIgY29uIHRyYW5zZm9ybWFjacOzbiBsb2dhcml0bWljYQoKI01vZGVsbyBNQSgxKSBlc3RhY2lvbmFsCm1vZGVsbyA9IEFyaW1hKEFpclBhc3NlbmdlcnMsIGMoMCwgMSwgMSksc2Vhc29uYWwgPSBsaXN0KG9yZGVyID0gYygwLCAxLCAxKSwgcGVyaW9kID0gMTIpLGxhbWJkYSA9IDApCmNvZWZ0ZXN0KG1vZGVsbykKbW9kZWxvYWx0ZXI9IEFyaW1hKEFpclBhc3NlbmdlcnMsIGMoMSwgMSwgMCksc2Vhc29uYWwgPSBsaXN0KG9yZGVyID0gYygxLCAxLCAwKSwgcGVyaW9kID0gMTIpLGxhbWJkYSA9IDApCgojIyBBbsOhbGlzaXMgZGUgcmVzaWR1YWxlcwojeDExKCkKcmVzaWR1YWxlcyA8LSBtb2RlbG8kcmVzaWR1YWxzCnBsb3QocmVzaWR1YWxlcykKYWNmKHJlc2lkdWFsZXMsbGFnLm1heCA9IDI0KQpwYWNmKHJlc2lkdWFsZXMsbGFnLm1heCA9IDI0KQojVGVzdCBkZSBhdXRvY29ycmVsYWNpP24KQm94LnRlc3QocmVzaWR1YWxlcywgbGFnID0gKGxlbmd0aChyZXNpZHVhbGVzKS80KSwgdHlwZSA9ICJManVuZy1Cb3giLCBmaXRkZiA9IDIpCiMjIyMjI0Fuw6FsaXNpcyBkZSBPdXRsaWVycwojVGVzdCBkZSBub3JtYWxpZGFkCmphcnF1ZS5iZXJhLnRlc3QocmVzaWR1YWxlcykKCgoKIyMjRXN0YWQ/dGljYXMgQ1VTVU0KcmVzPXJlc2lkdWFsZXMKY3VtPWN1bXN1bShyZXMpL3NkKHJlcykKTj1sZW5ndGgocmVzKQpjdW1xPWN1bXN1bShyZXNeMikvc3VtKHJlc14yKQpBZj0wLjk0OCAjIyNDdWFudGlsIGRlbCA5NSUgcGFyYSBsYSBlc3RhZD9zdGljYSBjdXN1bQpjbz0wLjE0NDIyIyMjI1ZhbG9yIGRlbCBjdWFudGlsIGFwcm94aW1hZG8gcGFyYSBjdXN1bXNxIHBhcmEgbi8yCkxTPUFmKnNxcnQoTikrMipBZipjKDE6bGVuZ3RoKHJlcykpL3NxcnQoTikKTEk9LUxTCkxRUz1jbysoMTpsZW5ndGgocmVzKSkvTgpMUUk9LWNvKygxOmxlbmd0aChyZXMpKS9OCnBhcihtZnJvdz1jKDIsMSkpCnBsb3QoY3VtLHR5cGU9ImwiLHlsaW09YyhtaW4oTEkpLG1heChMUykpLHhsYWI9InQiLHlsYWI9IiIsbWFpbj0iQ1VTVU0iKQpsaW5lcyhMUyx0eXBlPSJTIixjb2w9InJlZCIpCmxpbmVzKExJLHR5cGU9IlMiLGNvbD0icmVkIikKI0NVU1VNIFNxdWFyZQpwbG90KGN1bXEsdHlwZT0ibCIseGxhYj0idCIseWxhYj0iIixtYWluPSJDVVNVTVNRIikgICAgICAgICAgICAgICAgICAgICAgCmxpbmVzKExRUyx0eXBlPSJTIixjb2w9InJlZCIpICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgCmxpbmVzKExRSSx0eXBlPSJTIixjb2w9InJlZCIpCgpgYGAKCiMjIFByb27Ds3N0aWNvCgpgYGB7ciBQcm9ub3N0aWNvfQojeDExKCkKUHJvbm9zdGljb3M9Zm9yZWNhc3QobW9kZWxvLGg9MTIsbGV2ZWw9MC45NSkKcGxvdChQcm9ub3N0aWNvcykKcHJlZGljPC1wcmVkaWN0KG1vZGVsbyxuLmFoZWFkPTEyKQpwbG90KHByZWRpYyRwcmVkKQoKCiMjIyMjQ29tcGFyYWNpw7NuIGRlIHByb27Ds3N0aWNvcyMjIyMKbGlicmFyeShmcHApCnRyYWluIDwtIHdpbmRvdyhBaXJQYXNzZW5nZXJzLHN0YXJ0PWMoMTk0OSwwMSksZW5kPWMoMTk1OSwxMikpCnRlc3QgPC0gd2luZG93KEFpclBhc3NlbmdlcnMsc3RhcnQ9YygxOTYwLDAxKSxlbmQ9YygxOTYwLDEyKSkKZml0bW9kZWxvIDwtIEFyaW1hKEFpclBhc3NlbmdlcnMsIGMoMCwgMSwgMSksc2Vhc29uYWwgPSBsaXN0KG9yZGVyID0gYygwLCAxLCAxKSwgcGVyaW9kID0gMTIpLGxhbWJkYSA9IDApCnJlZml0IDwtIEFyaW1hKEFpclBhc3NlbmdlcnMsIG1vZGVsPWZpdG1vZGVsbykKZmMgPC0gd2luZG93KGZpdHRlZChyZWZpdCksIHN0YXJ0PWMoMTk2MCwxKSkKCgpoIDwtIDEKdHJhaW4gPC0gd2luZG93KEFpclBhc3NlbmdlcnMsc3RhcnQ9YygxOTQ5LDAxKSxlbmQ9YygxOTU5LDEyKSkKdGVzdCA8LSB3aW5kb3coQWlyUGFzc2VuZ2VycyxzdGFydD1jKDE5NjAsMDEpLGVuZD1jKDE5NjAsMTIpKQpuIDwtIGxlbmd0aCh0ZXN0KSAtIGggKyAxCmZpdG1vZGVsbyA8LSBBcmltYShBaXJQYXNzZW5nZXJzLCBjKDAsIDEsIDEpLHNlYXNvbmFsID0gbGlzdChvcmRlciA9IGMoMCwgMSwgMSksIHBlcmlvZCA9IDEyKSxsYW1iZGEgPSAwKQpmYyA8LSB0cyhudW1lcmljKG4pLCBzdGFydD1jKDE5NjAsMDEpLCBmcmVxPTEyKQpmb3IoaSBpbiAxOm4pCnsgIAogIHggPC0gd2luZG93KEFpclBhc3NlbmdlcnMsIGVuZD1jKDE5NTksIDEyKyhpLTEpKSkKICByZWZpdCA8LSBBcmltYSh4LCBtb2RlbD1maXRtb2RlbG8pCiAgZmNbaV0gPC0gZm9yZWNhc3QocmVmaXQsIGg9aCkkbWVhbltoXQp9CmRpZmU9KHRlc3QtZmMpXjIKZWNtPSgxLyhsZW5ndGgodGVzdCkpKSpzdW0oZGlmZSkKZWNtCmBgYAoKYGBge3IgT3Ryb3MgTW9kZWxhbWllbnRvfQpsaWJyYXJ5KHRpZHl2ZXJzZSkKbGlicmFyeSh0aWR5cXVhbnQpCmxpYnJhcnkoVFNzdHVkaW8pCmxpYnJhcnkoU0xCREQpCmxpYnJhcnkodGltZXRrKQpsaWJyYXJ5KGZhYmxlKQpsaWJyYXJ5KGZlYXN0cykKQWlyUGFzc2VuZ2Vyc190c2I9YXNfdHNpYmJsZShBaXJQYXNzZW5nZXJzKQpBanVzdGVfdHNpYmJsZV9QYXNzPC0gQWlyUGFzc2VuZ2Vyc190c2IgJT4lCiAgbW9kZWwoCiAgICBzYXJpbWEwMTEwMTFfUGFzcyA9IEFSSU1BKGxvZyh2YWx1ZSkgfiBwZHEoMCwxLDEpICsgUERRKDAsMSwxKSksCiAgICBhdXRvX3Bhc3MgPSBBUklNQShsb2codmFsdWUpLCBzdGVwd2lzZSA9IFRSVUUsIGFwcHJveCA9IEZBTFNFKSwgCiAgZHVtbXlfUGFzcz1BUklNQShsb2codmFsdWUpIH4gcGRxKDAsMSwxKSsgc2Vhc29uKCkpLAogZm91cmllcl9QYXNzPUFSSU1BKGxvZyh2YWx1ZSkgfiBwZHEoMCwxLDEpKyBmb3VyaWVyKEsgPSAyKSkKICkKCmdsYW5jZShBanVzdGVfdHNpYmJsZV9QYXNzKQoKI0FqdXN0ZV90c2liYmxlX1Bhc3MgJT4lIHNlbGVjdChzYXJpbWEwMTEwMTFfUGFzcyklPiUKICMgZ2dfdHNyZXNpZHVhbHMoKQoKI2F1Z21lbnQoQWp1c3RlX3RzaWJibGVfUGFzcykgJT4lCiMgIGZpbHRlcigubW9kZWw9PSdzYXJpbWEwMTEwMTFfUGFzcycpICU+JQojICBmZWF0dXJlcyguaW5ub3YsIGxqdW5nX2JveCwgbGFnID0gMTAsIGRvZiA9IDMpCgpBanVzdGVfdHNpYmJsZV9QYXNzICU+JQogIGZhYmxldG9vbHM6OmZvcmVjYXN0KGg9MTIpICU+JQogIGF1dG9wbG90KEFpclBhc3NlbmdlcnNfdHNiKQoKQWp1c3RlX3RzaWJibGVfUGFzcyAlPiUKICBmYWJsZXRvb2xzOjpmb3JlY2FzdChoPTEyKSAlPiUKICBmaWx0ZXIoLm1vZGVsPT0nc2FyaW1hMDExMDExX1Bhc3MnKSAlPiUKICBhdXRvcGxvdChBaXJQYXNzZW5nZXJzX3RzYikKCkFqdXN0ZV90c2liYmxlX1Bhc3MgJT4lCiAgZmFibGV0b29sczo6Zm9yZWNhc3QoaD0xMikgJT4lCiAgZmlsdGVyKC5tb2RlbD09J2R1bW15X1Bhc3MnKSAlPiUKICBhdXRvcGxvdChBaXJQYXNzZW5nZXJzX3RzYikKCkFqdXN0ZV90c2liYmxlX1Bhc3MgJT4lCiAgZmFibGV0b29sczo6Zm9yZWNhc3QoaD0xMikgJT4lCiAgZmlsdGVyKC5tb2RlbD09J2ZvdXJpZXJfUGFzcycpICU+JQogIGF1dG9wbG90KEFpclBhc3NlbmdlcnNfdHNiKQoKYGBgClZhbGUgbGEgcGVuYSBkZWNpciBxdWUgc2kgbG9zIHJlc2lkdWFsZXMgc29uIG5vIGNvcnJlbGFjaW9uYWRvcyBwZXJvIGVsbG9zIG5vIHRpZW5lIGRpc3RyaWJ1Y2nDs24gbm9ybWFsLCBzZSBwdWVkZSB1c2FyIGxvcyBpbnRlcnZhbG9zIGJvb3RhdHJhcCB2ZXIgc2VjY2nDs24gNS41IGxpYnJvIGZwcDMuIEVzdG8gc2UgaGFjZSBhw7FhZGllbmRvIGJvb3RzdHJhcD1UUlVFIGVuIGxhIGZ1bmNpw7NuICBmb3JlY2FzdCgpLgoK